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SINGLE -THERMOCOUPLE METHOD FOR DETERMINING 


HEAT FLUX TO A THERMALLY THICK WALL 

By Floyd G. Howard 
Langley Research Center 

SUMMARY 

A method that accounts for variable thermal properties has been developed for 
determining heat flux from a single temperature-time history within a thermally thick 
wall. The accuracy depends on the thermal depth of the measurement point, and ranges 
from good to excellent for a number of practical cases, with best results occurring for 
temperature measurements nearest the heated surface in conjunction with a small com- 
puting interval. The method is believed to be as accurate, if not more so, than other 
available methods utilizing multiple thermocouples. If measurement is near the heated 
surface and accurate, the method can detect sudden changes in slope of heat flux such as 
are experienced when boundary -layer flow becomes transitional or when heating-rate 
history fluctuates as a result of body motions in flight. 

INTRODUCTION 

One problem of considerable importance in the field of heat transfer is that of 
determining, with a reasonably high degree of accuracy, the heat flux to a calorimeter 
from measured temperatures. This problem is commonly referred to as the inverse 
problem. When the total (time-integrated) heat load permits the use of a thin-wall calo- 
rimeter, the task of determining the heat -flux history is not a very difficult one. How- 
ever, when test conditions necessitate a thick-wall calorimeter, the task becomes much 
more complex. In the thick-wall calorimeter, the variations through the wall of the 
time derivative of the temperature, the temperature itself, and the temperature-dependent 
properties (namely conductivity and specific heat) are significant and cannot be neglected 
as in the thin-wall calorimeter. 

A number of methods for calculating heating rates from measured temperatures 
have been developed (see refs. 1 to 12). However, most of these solutions assume con- 
stant thermal properties and therefore have inherent inaccuracies in most practical 
cases, where temperatures vary over a wide range. The method developed in refer- 
ence 2 for the reduction of Project Fire calorimeter data does account for variable ther- 
mal properties; however, it does have some disadvantages, the major one being that it 



requires data from several thermocouples in depth. The disadvantage of this method, 
when a practical number of thermocouples (like four) are used, is that all thermocouples 
must function properly in order to obtain accurate data. 

The realization of the need for an alternate method that could render valid data 
analysis if temperature data were not obtained from all thermocouples prompted an inves- 
tigation which led to the development of the new method reported herein. This method 
can determine heat flux from the temperature -time history measured at only one depth. 
Therefore, data from a greater number of stations can be obtained when the total number 
of thermocouples is limited. The method accounts for variable thermal properties (as 
functions of temperature) and can account for the effects of radiation losses or a com- 
posite wall. It is believed to be as accurate if not more so than other known methods 
and can probably be modified to account for a coolant gas being blown through a porous 
wall. 


SYMBOLS 

A area, centimeter^ 

A square coefficient matrix 

B column vector 

c specific heat of wall material, joules/gram -°K 

E t temperature error (T com - T input ), °K 

k thermal conductivity of wall material, watts/centimeter-°K 

q heat flux, watts/centimeter^ 

q^ erroneous heat flux, watts/centimeter^ 

R radius of Fire II forebody calorimeter shield (ref. 3) measured perpendicularly 

from axis of symmetry, centimeters 

s radial distance measured perpendicularly from axis of symmetry of Fire II 

reentry vehicle (ref. 3), centimeters 

T temperature, °K 
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T’ temperature at beginning of computing time interval, °K 

t time, seconds 

At computing time interval, seconds 

X unknown column vector 

e emissivity 

p density of wall material, grams/centimeter ^ 

a degree of data scatter; Stefan-Boltzmann constant, 

5.6697 x 10~12 watts/centimeter2-(°K)^ 

r depth from heated surface to center of block, centimeters 

At block thickness in thermal model, centimeters 

r total wall thickness, centimeters 

Subscripts: 

1,2, 3, 4 thermocouple (TC) number; time-interval number; block number 
com computed on basis of assumed q 

input input or desired 

m mean value 

PROCEDURE AND MECHANICS OF NEW METHOD 

Solution for the heat flux is found by using a direct thick-wall heating computer 
program which computes temperatures for a given heating-rate history (see appendix). 
The wall thickness is assumed to be divided into a number of blocks, and a finite- 
difference computational technique is used to compute the temperature of each block. 

A modification to this program provides for iteration on an assumed heating rate during 
a short computing time interval to find the value required to produce a known temperature 
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rise at a particular depth. The iteration is very rapid since the computed temperature 
error (i.e., the difference between the computed temperature rise based on the assumed 
heat flux and the desired temperature rise) will vary linearly with the magnitude of the 
assumed heat flux. Therefore, a linear interpolation or extrapolation of two assumed 
values of heat input and their associated temperature errors provides the required value 
of heat input for temperature error of zero, as illustrated in detail in figure 1. 

The solid line in the plot of temperature as a function of time (fig. 1(a)) represents 
a known temperature-time history at some known depth in the wall (r). First, an assumed 
heat flux (indicated by the circle marked 1 in fig. 1(b)) is used to compute the tempera- 
ture rise over the time interval. This temperature rise (indicated by the circle marked 1 
in fig. 1(a)) is in error Erp at the end of the time interval. This temperature error as 
a function of the assumed heat flux is indicated by the circle marked 1 in figure 1(c). 

The procedure is repeated with a second assumed value of heat flux (indicated by the 
square marked 2 in fig. 1(b)) to provide a second temperature error as a function of heat 
flux. This second temperature error is indicated by the square marked 2 in figure 1(c). 

If thermal properties are held constant for the small computing time interval (At) as is 
done in the finite -difference heating program, the heat flux varies linearly with the tem- 
perature error at any depth in the wall. Therefore, two assumed values of heat flux and 
their associated temperature errors can provide a third and correct value of q at 
= 0 by linear interpolation or extrapolation. In order to provide the correct tem- 
perature for each block in the wall, the computations of temperature rise are repeated 
with this correct heat flux used. The correct temperature is then used to start the com- 
putations for the next time interval. This procedure is repeated in order to obtain a 
solution for q for each time interval. In the example shown in figure 1, the solutions 
determined from the heat -flux assumptions 1 and 2, 4 and 5, and 7 and 8 are indicated by 
the solid triangles marked 3, 6, and 9, respectively, in figures 1(b) and 1(c). Since the 
heating rates are held constant over the computing time interval in the finite-difference 
computations in the thick -wall program, the solution obtained is a step-changing function 
as indicated by the heavy, solid horizontal lines in figure 1(b). In order to obtain a 
smoothly changing value of heat rate (q), the solution for each time interval should be 
considered to occur as an instantaneous value at the center of the interval, as indicated 
by the line faired through the solid triangles in figure 1(b). With the exception of the 
first time interval, the first assumed value of q is always taken to be a 25-percent 
increase over the correct value determined for the previous time interval (the user must 
provide the first assumption for the very first interval in order to start the computations). 
The second value of q is always a 2 5 -percent increase over the first assumption. This 
increase of 25 percent, however, is arbitrary. Changing this arbitrary increase to 5 per- 
cent was investigated and found to have no effect. 
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The basic inputs required for the single -thermocouple inverse program, in addition 
to those required for the direct thick-wall heating program (initial temperature gradient 
through the wall, material properties, dimensions, etc.), are the measured temperature- 
time history, the depth of the measurement (identified by the number of the block that has 
the proper depth in the thermal model), and the initial assumed heat flux for the first com- 
puting interval. The outputs are the computed heat -flux history and the temperature his- 
tory for each block. 

If the thick wall of a flight vehicle is subjected to heating before the start of compu- 
tations, the initial temperature gradient through the wall will not be defined by only one 
thermocouple in depth. If it is desirable to start the inverse computations at some time 

when ( ^~ ^ 0 (and therefore ^ 0 ], a temperature gradient through the wall can be 
dt \ dr / 

assumed. The effects of an erroneous but reasonable initial temperature gradient will 
diminish very quickly. 

All the computations shown in the present report were obtained by the technique 
described in this section. However, based on the findings in this report, a more recent 
program that solves directly for the heat flux without an iteration has been written. The 
reason q can be solved for directly is explained briefly in the appendix. This factor, 
however, is important only in regard to machine time and has no significant bearing on 
the accuracy of the computations. 

ANALYSIS AND RESULTS 
Analytical Investigation 

The validity and accuracy of the single-thermocouple inverse method were analyti- 
cally investigated by comparing results from the method with known heat inputs. This 
investigation was conducted by computing the temperatures through a 1.524-cm-thick 
beryllium wall, divided into 30 blocks to provide an accurate thermal model. The ther- 
mal properties used for beryllium are shown in figure 2. Figure 3 shows the geometrical 
details of the 1.524-cm-thick thermal model, namely the thickness of each block and the 
depth of the block center or node from the heated surface. Thinner blocks are used near 
the heated surface because of the large temperature gradients which exist in this region 
at high heating values. In order to obtain the temperature close to the back surface, thin 
blocks were used in this region also. With an assumed heat -flux history as input, the 
basic direct thick-wall heating computer program provided a temperature history for 
each block. Then it was assumed that a temperature -history measurement was made at 
only one depth in the wall. This temperature history and the depth at which it occurred, 
specified by block number, were used as inputs to the single-thermocouple inverse pro- 
gram. The heat flux computed by the inverse method was compared with the input heat 
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flux. This comparison was made for three different types of heat -flux histories — 
namely, a smoothly changing one representative of a reentry flight test at zero angle of 
attack, a step-input heat -flux history representative of a test in a ground facility, and an 
oscillating heat -flux history that might result from a flight test with significant body 
motions. 

In order to ascertain that the thermal model used in the analytical investigation 
was divided into a sufficient number of blocks (30), the effect of the number of blocks on 
accuracy was investigated. 

First, the smoothly changing heat -flux history was used in the direct program with 
an 11-block model, a 30-block model, and a 70-block model. All models were 1.524 cm 
thick, and the computing time interval was 0.05 second. Temperature-time histories at 
identical depths were compared and found to be in very close agreement. For example, 
at 13 seconds the temperature at the 0.254-cm depth was 1424.77° Kfor the 11-block 
model and 1420.67° K for the 30-block model after a rise of over 1080° K. Because of 
the long computer time required for the 70 -block model, the run for that model was ter- 
minated at 4.4 seconds. The temperature at the 0.254-cm depth was 387.34° K, 387.46° K, 
and 387.56° Kfor the 70-, 30-, and 11-block models, respectively, after a rise of over 
53° K. 

Next, in order to be sure that a 30 -block model was also sufficient for the single - 
thermocouple inverse solution, the temperature -time history at the 0.254-cm depth was 
used in inverse solutions for the 11 -block thermal model and then for the 30-block model. 
The computed heat -flux histories obtained by using the two different thermal models and 
a computing time interval of 0.1 second in the inverse solution (with the use of the same 
temperature history) differed by less than 0.3 percent. It was concluded that the 30-block 
model gave sufficiently fine thermal detail and that any increase in the number of blocks 
would have insignificant effect on the results. 

The first of the three different heat -flux histories used in the investigation is shown 
in figure 4. This history represents the heating on a reentry-flight -test vehicle with zero 
trim and no oscillatory angle of attack. The time at which boundary-layer transition was 
assumed to begin is noted. The temperatures at four depths that result from this heating 
are also shown in the figure. 

It should be noted that in figure 5, which shows solutions by the single -thermocouple 
inverse method, beryllium would start to melt at the surface at approximately 11.8 sec- 
onds for the heat flux shown. However, in order to analyze the method during and after 
the time of maximum heat flux, a hypothetical case has been assumed in which melting 
has been disregarded and thermal properties above melting temperature have been held 
constant (held at melting -temperature values). 


6 



Figure 5(a) shows the effect of computing interval (0.1, 0.2, 0.25, and 0.5 second) 
at a measurement depth r of 0.0254 cm. The agreement between the computations and 
the exact input becomes better as the computing interval is decreased from the largest 
interval to the smallest interval. Figure 5(b) shows the effect of computing interval 
(0.2 and 0.5 second) at a measurement depth t of 0.2540 cm. At this deeper location, 
computations by use of the 0.2-second interval show a tendency to oscillate during 
increasing heat flux. The oscillations become more noticeable after maximum heating. 
Although they are not shown in figure 5(b), computations for a computing interval of 
0.1 second at t = 0.2540 cm were made and found to oscillate a little more than those 
for a computing interval of 0.2 second, whereas computations for the 0.5-second interval 
show no indication of oscillation and are in good agreement with the exact input. In the 
finite-difference computer program the thermal properties (which actually vary instanta- 
neously with temperature) are held constant over the computing interval at their values 
corresponding to temperature at the beginning of the interval (T’)« This practice obvi- 
ously is a source of error (which is small for examples shown). Therefore when temper- 
ature (and thus thermal properties) are changing rapidly, it is apparent that the accuracy 
will improve with the use of small computing intervals. However, for the finite- 
difference computing technique, oscillations in the computed heat flux result if the com- 
puting interval is too small for a particular thermal depth (dependent on measurement 
depth and thermal properties). 

The frequency at which the temperatures are provided also has influence on the 

magnitude and pattern of the oscillations because constant slope (— ] is assumed between 

\dt/ 

the provided points. Figure 6 shows how the computing interval and the frequency at 
which temperature is provided affect the oscillation of computations for a depth r of 
0.6096 cm (40-percent depth) in a 1.5 24 -cm -thick wall. Computations based on tempera- 
tures provided every 0.4 second are shown by the circles for computing intervals of 0.1, 
0.2, and 0.4 second in figures 6(a), 6(b), and 6(c), respectively. Oscillations occur for 
both the 0.1- and 0.2 -second computing intervals and tend to diverge at maximum heating. 
No significant oscillations occur for the 0.4-second computing interval (which corre- 
sponds to the frequency at which temperatures were provided). However, some disagree- 
ment (5 -percent maximum difference between the computed values and the exact input) is 
noticeable. (It should be noted that the heat -flux scale starts at 500 watts/cm^.) Results 
of computations with temperature inputs at 0.1- and 0.2 -second intervals, corresponding 
to the computing intervals, are shown by the square symbols in figures 6(a) and 6(b). The 
oscillations are reduced significantly for the computing interval of 0.1 second and are 
essentially eliminated up to time of maximum heat flux for the computing interval of 
0.2 second. Apparently, oscillations in solutions for q are sensitive to the magnitude 
dT 

of the steps of — at the times of temperature inputs (particularly for measurements 
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at large depths) and therefore can be minimized by providing temperature inputs for 
every computing interval. 



It should be realized that an oscillating heat flux to the surface may not cause 
oscillations in the temperature history at a depth in the wall (not at the heated surface). 
Therefore, it is believed that an absolute unique solution cannot be obtained without sur- 
face temperature, since an oscillating heat flux or a smoothly changing one may produce 
the same temperature -time history at a measurement depth not at the heated surface. 

Only at the heated surface, where the response time is zero, does a temperature history 
dictate one and only one heat -flux history, and many solutions which oscillate about some 
mean solution may exist for measurements at locations other than the heated surface. 

The oscillations tend to diminish or become less pronounced as the thermal depth becomes 
less. Therefore from this point of view, the most desirable thermocouple location is as 
near the heated surface as possible. Presently no suitable governing criteria have been 
established to determine the smallest computing time interval that can be used for a 
given thermal depth without causing oscillations in the solution. It may be noted that the 
lack of an absolute unique solution without surface measurement is not peculiar to the 
single-thermocouple inverse method. 


Shown in figure 7 are solutions for a constant or step-input heat flux of 
600 watts/cm^. A constant heat flux is characteristic of a ground facility. The initial 
assumed heat flux for the first computing interval was 200 watts/cm^. Solutions are for 
computing intervals of 0.2, 0.4, 0.6, and 0.8 second for each of three depths, 0.0254, 
0.2540, and 0.5080 cm. The accuracy is good even for a depth of 0.5080 cm except for a 
transient period, which is less than 0.5 second. For each of the three depths the accu- 
racy is shown to improve as the computing interval is decreased. One reason is that, 
as previously mentioned, the thermal -property values for each block are based on the 
temperature of the block at the beginning of the computing interval and are held constant 
over the interval. This procedure introduces a source of error, the magnitude of which 
depends on how much the block temperatures vary over the interval and how much the 
properties vary with temperature. It is apparent that this source of error could be 
diminished by extrapolating the block temperatures with time and basing the properties 
on the value of temperature at the middle of the computing interval. 

Because of body motions, flight vehicles may experience an oscillating heat flux 
rather than a smoothly changing one. The method was analyzed for such a case and the 
results are shown in figure 8. The oscillating heat-flux history used as the input was 
determined as follows: (1) the most severe body motions expected of a planned reentry 

vehicle were determined by a six-degree-of -freedom program; (2) from the body 
motions, an angle -of -attack history for a particular body location was determined; and 
(3) the effect of angle of attack on heating was determined by theoretical computations 


8 



and used to modify the zero-angle -of -attack heating (shown in fig. 5) and thus produce the 
oscillating heat-flux history shown in figure 8. Solutions by the single-thermocouple 
inverse method shown in figure 8 are for computing intervals of 0.05 and 0.10 second at 
a measurement depth of 0.0254 cm. At this depth there were noticeable although slight 
oscillations in the temperature history. At a depth of 0.254 cm the temperature history 
was even less responsive (as would be expected) to the oscillating heat flux; that is, the 
temperature histories became smoother as depth was increased. Therefore no attempt 
was made to reproduce the exact oscillating heat flux from a temperature at this depth. 
The oscillating temperatures at the 0.0254-cm depth were provided at 0.05-second inter- 
vals. The agreement shown in figure 8 by the comparison with the exact input indicates 
the ability of the method to detect sudden changes in the slope of the heat flux if the mea- 
surement is near enough to the heated surface that it will respond to the oscillations in 
the heat flux. It is very important to note that the temperatures put into the inverse 
solutions for figures 5 to 8 were precise (except for the linear variation assumed between 
points) — that is, without scatter or bias. This fact is especially important for the oscil- 
lating heat flux (fig. 8), since the temperature history at the 0.0254-cm depth had only 
slight oscillations. The oscillations were within a total bandwidth of less than 10° K 
during the period of maximum oscillation. It is believed that with current means of 
obtaining temperature measurements during flight, the oscillations may be undefinable 
because of scatter, for a case similar to the one shown. However, it should be understood 
that the shortcoming would be in defining the correct temperature history and not in the 
single -thermocouple inverse method. 

The results shown in figures 5, 7, and 8 indicate that the accuracy of the method 
ranges from good to excellent, with best results obtained for temperature measurements 
nearest the heated surface in conjunction with a small computing interval. Although 
these results were all obtained for a 1.524-cm-thick beryllium wall, it is believed that 
the trends would apply to other materials and wall thicknesses. No results that include 
effects of radiation losses or a composite wall are shown. However, the method can 
account for these effects. 


Comparison With Computation by Other Methods 

Com parison with integral method of reference 2 .- The rate at which heat is stored 
in a thick wall is equal to the density of the material (p) multiplied by the product of 

specific heat and change in temperature with respect to time ^c integrated from the 

heated surface to the back surface (ref. 13). If the density of the wall material is treated 
as a constant, the one-dimensional heat flux can be expressed as 
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9=p C( c f) di (1) 

The method developed in reference 2 approximates the exact solution of equation (1) 

dT 

by: (1) determining the value of c — from measured temperatures at four depths in 

dt 

the wall; (2) fitting an analytical expression through the four values; and (3) integrating 

dT 

analytically under the curve to obtain q. The four values of c — at a particular 

dt 

time (t) indicated by the circles in figure 9 are determined by fitting an analytical expres- 
sion f(t) to the measured temperature -time history, differentiating the expression to 
dT 

obtain — , and obtaining c as a function of the temperature. The accuracy in com- 
dt 

puting the heat flux by this method is dependent upon the accuracy of the analytical curve 

dT 

fit of both temperature as a function of time and c — as a function of t. This method 

dt 

of approximation by analytical integration (ref. 2) is hereinafter referred to as the inte- 
gral method. 

The squares in figure 10 represent computations by the integral method for a 
1.524-cm-thick beryllium wall (three times the thickness used in refs. 2 and 3) with mea- 
surement depths of 0.0254, 0.5080, 1.0160, and 1.5176 cm. At the lower heating rates the 
computations agree with the exact input very well but at the higher heating rates, differ- 
ence between the computed value and the exact input becomes noticeable. The measure- 
ment locations were spaced approximately evenly through the wall. This evenly spaced 
arrangement is the same as that in references 2 and 3. Because of the large gradients 

^variation of T and c see fig. llj near the heated surface, it was believed that 

with the 1.524-cm-thick wall, better results could be achieved by placing the second and 
third measurement locations nearer the heated surface. The diamonds in figure 10 rep- 
resent computations by the integral method for measurement depths of 0.0254, 0.3556, 
0.8636, and 1.5176 cm. At the lower heating rates the results are essentially the same 
as for the evenly spaced arrangement. At the higher heating rates the agreement between 
the computations and the input is somewhat worse than for the evenly spaced arrangement. 
The source of the error was investigated and found to be associated with the shape of the 
analytical curve (a third-order polynomial as in refs. 2 and 3) used to curve fit the four 

dT 

measured values of c — through the wall. As shown in figure 11, there is noticeable 

dt 

dT 

disagreement between the third-order polynomial and the exact variation of c — 

through the wall. In regard to the final answer, the important thing is, of course, not the 
accuracy in the shape of the polynomial but in the integrated area under the curve. For 
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the two sets of measurement depths shown, it is obvious that the shape of the polynomial 
more closely approximates that of the exact curve when measurement locations are 
evenly spaced. Although it is not so obvious, the integrated area for the evenly spaced 
arrangement is in a little better agreement than for the unevenly spaced arrangement. 
These results are apparent in the comparison of computed heat flux with exact input 
shown in figure 10. Although not shown in figure 10 or figure 11, results were deter- 
mined also for measurement depths of 0.0254, 0.2032, 0.6096, and 1.524 cm. The error 
in the area under the polynomial became magnified even more. At 11 seconds the com- 
puted heat flux was 1220 watts/cm^, which is 19.4 percent higher than the exact-input 
value of 1022 watts/cm^ shown in figure 10. Errors due to this source become smaller 
relative to the heat flux for smaller physical thickness of the wall. Therefore the influ- 
ence of the foregoing error on the results of references 2 and 3 was small, since the wall 
thickness was only 0.508 cm and the thermocouples were evenly spaced. 

Computations by the integral method shown in figure 10 start after transition and 
end before peak heating. The reason is that in the integral program (ref. 2) the 
temperature-time history was expressed as a sixth-order polynomial. Sixth-order poly- 
nomials fitted to the entire temperature history from 0 to 16 seconds for each of the four 
depths give very poor computed heat-flux results. Therefore, care must be taken to fit 
the temperature history with the polynomial over regions that have no drastic changes or 
inflection points as occur at transition and peak heating, respectively. It may be noticed 

dT 

in figure 11 that the values of c — - indicated by the symbols may not fall precisely on 

dT 

the exact curve. The reason is that the value of c — — determined from the sixth-order 

dt 

polynomial, which is fitted to the temperature history, may have a slight error. Although 

results for lower heat flux are not shown in figure 11, calculations at 6 seconds 

(q = 352 watts/cm2) were made which show that the agreement between the third-order 

dT 

polynomial and the exact curve of c — through the wall is very good (see agreement 

dt 

of q in fig. 10) when the gradient through the wall is small. Although a third-order 

dT 

polynomial fitted to c — as a function of r gave satisfactory results for the 0.508-cm- 

thick wall of references 2 and 3 (even for high heating levels), it is apparent from fig- 
ure 11 that some other analytical expression may give more accurate results for thick 
walls with high heating. Unpublished results of an investigation of this problem by 
James L. Dillon of Langley Research Center indicate that a third-order exponential pro- 
vides improved accuracy. 

Results determined by the single -thermocouple inverse solution for a depth of 
0.0254 cm and a computing interval of 0.1 second are also shown in figure 10 for 
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comparison. The agreement is good over the entire heating period. The input tempera- 
ture history was smooth (without scatter or bias) for both methods in figure 10. 

Comparison with sin gle -thermocouple method of reference 4.- The method of ref- 

dT 

erence 4 is based on the assumption that the mean value of c — always occurs at a 

dt 

specific depth in the wall, namely at approximately 40 percent of the total thickness mea- 
sured from the heated surface [ — = 0.40). It is obvious that there is a mean value of 

^ l T t I 

c — somewhere in the wall which if treated as a constant from the heated surface to 
dt 

the back surface (see dashed line in fig. 11) will give the same area as the integrated 

dT 

area under the exact curve. When c — - is treated as a constant, equation (1) is sim- 

dt 

plified and becomes 


Q = 



T t 


( 2 ) 


This method was used to determine heating rates from temperatures for the case 
representing a reentry flight experiment with wall thickness of 1.524 cm. Results com- 
puted by this method for the 40-percent depth (t = 0.6096 cm) are in very good agreement 
with the exact input at the lower heating values (as shown in fig. 12), but not at the higher 
heating values, or in the region of maximum heating. The nondimensional location (t/tj.) 

dT 

of the mean value of c is not fixed; for instance, at higher heating levels the mean 
dT 

value of c — progresses from the 40-percent depth toward the heated surface while 
dt 

the heat flux is increasing. The location is apparently dependent on material properties, 
heat flux, and wall thickness. As with the integral method, the accuracy of this method 
is also improved as the physical thickness of the wall becomes smaller, and it may be 
noted that good results were obtained for a 0.508-cm -thick beryllium wall, as shown in 
reference 4. 

Also shown in figure 12 for comparison is the present single -thermocouple inverse 
solution for a measurement which is at the 40-percent depth. The computing interval is 
0.4 second. The results are significantly better than those obtained by the method of 
reference 4 (eq. (2)). 


Influence of Erroneous Inputs 

There are several sources of error in the determination of q (some of which are 
obvious), namely incorrect material properties, error in location of thermocouple, 
biased temperature history, and scattered temperature data. The effects of these factors 
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were investigated for a thermocouple location 0.0254 cm below the heated surface of a 
1.524-cm-thick beryllium wall (see block 2 in fig. 3). It should be noted that the errors 
shown in tables I to IV are due specifically to error in the parameter being considered; 
that is, they reflect only the effect of the parameter error and not the accuracy of the 
single -thermocouple inverse method. 

Error in specific heat .- In order to determine the effect of error in specific heat, 
an assumed error in the curve for c shown in figure 2 was used in the single- 
thermocouple inverse solution. The effects of errors in specific heat of -2.5 percent 
and -5 percent are shown in table I. The percent error in the heat flux is not the same 
as the percent error in the specific heat, as might be expected from equation (1). Appar- 
ently, the reason is that the error in the specific heat alters the change in temperature 

with respect to time through the wall. In this case it appears that the secondary 
dT 

effect on — tends to partially compensate the error in specific heat. The effect of 
dt 

this type error does not diminish as physical thickness of the wall becomes smaller. 

Error in conductivity .- In order to determine the effect of error in conductivity, 
an assumed error in the curve for k shown in figure 2 was used in the single - 
thermocouple inverse solution. The effect of errors in conductivity of -2.5 percent, 

-5 percent, -10 percent, and +10 percent is shown in table n. The effect on the computed 
heat flux is as expected. However, it should be noted that with the same error in conduc- 
tivity but with the measurement near the back wall, the effect on the computed heat flux 
would have been just the opposite in sign and that therefore there is a cross-over point 
somewhere within the wall where the effect of error in conductivity would be minimized. 
The effect of this type error diminishes as physical thickness of the wall becomes 
smaller. 


It should be noted that when two or more thermocouple measurements are used in 
combination (as discussed subsequently), the effect of erroneous conductivity tends to be 
minimized if the thermocouple distribution nearly spans the wall thickness. 


Error in measurement location.- Since the temperature and the temperature rise 


with respect to time 



vary with distance from the heated surface, it is essential that 


the location of the measurement be known with reasonable accuracy. Effects of location 
errors of 0.0254 cm and 0.0508 cm (too deep) are shown in table HI; thus, temperature 
histories at 0.0508-cm and 0.762-cm depths were considered to be actually at a depth of 
0.0254 cm. The error in the computed heat flux due to the larger location error of 
0.0508 cm (see results in table for temperature history at 0.0762 cm) is approximately 
double that due to the smaller location error of 0.0254 cm (see results for temperature 
history at 0.0508 cm). The computed heat-flux errors are comparatively small even for 
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these large location errors (locations can be determined in practice within about 
±0.008 cm). The effect of this type error diminishes as the physical thickness becomes 
smaller. 

Biased tempe rature.- Effect of temperature which is 5 percent too high is shown 
in table IV. It may be noted that most of the effect shown is due to the increase in rate 

of temperature rise associated with temperature that is too high by a constant per- 
centage. Temperatures too high by a constant number of degrees would have had a rela- 
tively small effect on the results. 

Data scatter .- All the aforementioned solutions were determined by using the 
smooth temperature-time histories computed by the direct thick -wall heating computer 
program. A solution was also obtained wherein random scatter, simulating random 
measurement and telemetry inaccuracies (no bias inaccuracies), was applied to the ini- 
tially smooth temperature -time history. A random-number generator was used to pro- 
vide statistical scatter for a 3a value of 25° K. The scattered temperature data were 
smoothed by a Chebyshev curve fit of a sixth-order polynomial which was then used in 
the single -thermocouple inverse program to obtain a heat -flux history. The heating 
results were found to be in very good agreement with the exact input (within 2.5 percent), 
as shown in figure 13. The sixth-order polynomial curve was fitted to the temperature 
history over the time range of 3 to 12 seconds, since fitting the curve from a time before 
the abrupt change at transition to a time after maximum heating proved to be a very poor 
representation of the real temperature -time curve. It may be noted that this computation 
is actually a test of the accuracy of the sixth -order -polynomial curve-fit technique rather 
than a test of the single -thermocouple method per se. 

SOME APPLICATIONS 
Application to Flight Data 

Heating-rate — time histories have been computed by the single -thermocouple 
inverse method for a location on each of the three beryllium heat shields (three data 
periods) of the Fire II reentry spacecraft (ref. 3). The thickness of each of the heat 
shields was 0.508 cm. The identical sixth -order polynomial curves which were fitted to 
the experimental temperature data in the Project Fire data analysis were used herein. 

As shown in figures 14 and 15, the stagnation-point heating rates during the first 
and third data periods, as computed by the single -thermocouple inverse method with 
either of the two outermost thermocouples, are in good agreement with the data from 
reference 3, determined by the integral method. Since the Fire n flight data are confi- 
dential, the scale values are not given in figures 14 to 18 in order that these data can be 
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included herein. Figure 16 shows results from the second data period for locations on 

the heat shield of 0.54 and 0.90 (the comparison was made at these locations 

rather than at the stagnation point because of failure of thermocouples at the stagnation 
point). Because of the noticeable disagreement between results from the two outermost 
thermocouples at both s/R locations, single-thermocouple-method results were com- 
puted by using each of the four thermocouples at these two stations. The integral- 
method data of reference 3 is roughly the average of the results from the four individual 
thermocouples. The source of the difference between results from the individual thermo- 
couples during this data period has not been determined, although certain possibilities 
have been considered. Temperature gradients through the wall at the beginning of the 
heating period, which were initially suspected as possible contributing factors, were 
studied and found to damp out quickly after the heating commenced. Also, the flight tem- 
perature measurements at different radial locations indicate that lateral conduction was 
negligible and could not be the source of the discrepancies. It may be noted, however, 
that the temperature rise was much more rapid during the second data period than during 
the first or third period, as shown in figure 17. This rapid rise results in fewer points 
per unit temperature rise and thus less accuracy in defining the temperature curves. 

For this reason more variation between the individual-thermocouple results might be 
expected during the second data period than during the other data periods. 

It may be noted that a machine technique has been developed to utilize multiple - 
thermocouple data with the basic, single -thermocouple inverse method. The same thick- 
wall, finite -difference computation program is used, but instead of finding the heating 
rate q that produces the one measured temperature -time history, the criterion for the 
desired q is that the sum of the temperature errors at the four measurement depths 
(or as many as are available) be zero. That is, the sum of the computed temperatures 
at the measurement depths must equal the sum of the experimental temperatures. This 
method is believed to be an improvement over averaging the values of q obtained from 
individual measurements containing experimental inaccuracies, since it effectively gives 
more weight to the more sensitive measurements (i.e., those thermocouples nearer the 
heated surface). The heating rates computed by this method in the second data period 

for the location —= 0.90 are compared in figure 18 with the integral-method results of 
R 

reference 3. The two sets of results are within 2 or 3 percent of each other. 

Other Possible Applications 

It is believed that the single -thermocouple inverse method could be used to account 
for gas (as a coolant) being blown through a porous wall. This application would require 
that the heat-balance equations (shown in the appendix) be extended to include the effects 
of the mass -flow rate and physical properties of the gas. 
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The method could also conceivably be used in other applications, such as deter- 
mining the conductivity (or specific heat) of a material by iterating, for example, an 
assumed value of conductivity until temperature measurements at two or more different 
depths are satisfied. 


CONCLUDING REMARKS 

A method for determining heating rates from a single temperature -time history 
within a thermally thick wall has been developed. This method accounts for variable 
thermal properties, can be used for a composite wall, and can account for surface- 
radiation effects in the analytical computations if so desired. The method has been 
analyzed for a 1.524 -cm -thick beryllium wall and the following points noted: 

1. The accuracy depends primarily on the thermal depth of the measurement point, 
and ranges from good to excellent for a number of practical cases, with best results 
occurring for temperature measurements nearest the heated surface in conjunction with 
a small computing interval. 

2. If the measurement is sufficiently near the heated surface, the method can detect 
sudden changes in the slope of heat flux such as are experienced when boundary -layer 
flow becomes transitional or when heating-rate history oscillates. 

3. A comparison with an integral method and with another single -thermocouple 
method shows that the present method is more accurate for the heat -flux histories and 
the thermal model considered herein. 

4. The method is not restricted to data from only one measurement, but can utilize 
any number of measurements and therefore all available data for any given station. 

Langley Research Center, 

National Aeronautics and Space Administration, 

Langley Station, Hampton, Va., May 14, 1968, 

711-02-09-01-23. 
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APPENDIX 


THICK-WALL HEAT -BALANCE EQUATIONS 

The basic thick-wall heating computer program utilizes a backwards finite- 
difference computing technique. The finite -difference equations are solved simulta- 
neously by the Gauss-Seidel method. The heat -balance equations used in the program 
for one -dimensional heat flow are (for a wall divided into four blocks) as follows: 


Heat Heat 

in radiated 

i out 



Heat stored 


Heat 

conducted out 


AATjpc 




(Al) 


Heat 

conducted in 


2A 




Heat stored 


AAT2PC 




Heat 

conducted out 
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Heat 
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Heat stored 
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Heat 

conducted out 



(A3) 


Heat 

conducted in 



Heat stored 
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Heat 

radiated 

out 
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APPENDIX 


In these equations T’ is the temperature at the beginning of the time interval and 
T is the temperature at the end of the time interval. Thermal properties (c, k, and e) 
are based on T'. 

In the basic thick -wall program operation, q is known and the unknowns are T^, 
Tg, Tg, and T^. Therefore, a unique solution exists for the unknown temperatures 
since the number of independent equations provided and the number of unknowns are the 
same. In the single -thermocouple inverse method as described in the section "Procedure 
and Mechanics of New Method," an assumed value of q is used in order to solve for the 
temperature at each block. If, however, the temperature history is known at one block, 
it is not necessary to solve for temperature at that location, and thus q could be con- 
sidered as one of the unknowns without increasing the number of unknowns over the num- 
ber of equations. Then q could be solved for directly without going through an iteration 
on an assumed q. A program utilizing a matrix solution, which will solve for q 
directly, has been written. This system is solved by a library subroutine which solves 
the matrix equation AX = B (where A is a square coefficient matrix and B is a 
column vector). The heating rate is evaluated by a subroutine in place of the known 
temperature. 
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TABLE I.- EFFECT OF ERROR IN SPECIFIC HEAT (c) ON 
INVERSE HEAT -FLUX COMPUTATIONS 


Time, 

sec 

<l2’ 

watts/cm2 

<lE’ 

watts/cm^ 

<lE “ 12’ 
watts/cm2 

Percent 

error 


Error in specific heat of - 

2.5 percent 


1.1 

26.39 

26.23 

-0.16 

-0.60 

2.1 

29.59 

29.17 

-.42 

-1.42 

3.1 

65.94 

64.90 

-1.04 

-1.57 

4.1 

170.19 

167.67 

-2.52 

-1.48 

5.1 

276.17 

271.93 

-4.24 

-1.53 

6.1 

363.92 

358.06 

-5.86 

-1.61 

7.1 

471.32 

463.56 

-7.76 

-1.65 

8.1 

595.30 

585.44 

-9.86 

-1.66 

9.1 

743.02 

730.77 

-12.25 

-1.65 

10.1 

891.44 

876.89 

-14.55 

-1.63 

11.1 

1043.89 

1027.07 

-16.82 

-1.61 


Error in 

specific heat of - 

■5 percent 


1.1 

26.39 

26.06 

-0.33 

-1.25 

2.1 

29.59 

28.74 

-.85 

-2.88 

3.1 

65.94 

63.86 

-2.08 

-3.15 

4.1 

170.19 

165.11 

-5.08 

-2.99 

5.1 

276.17 

267.64 

-8.53 

-3.09 

6.1 

363.92 

352.11 

-11.81 

-3.24 

7.1 

471.32 

455.69 

-15.63 

-3.32 

8.1 

595.30 

575.43 

-19.87 

-3.34 

9.1 

743.02 

718.36 

-24.66 

-3.32 

10.1 

891.44 

862.15 

-29.29 

-3.29 

11.1 

1043.89 

1010.05 

-33.84 

-3.24 
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TABLE H.- EFFECT OF ERROR IN CONDUCTIVITY (k) ON 
INVERSE HEAT -FLUX COMPUTATIONS 


Time, 

sec 

<l2’ 

watts/cm2 

<lE’ 

watts/cm2 

^E ~ q 2’ 
watts/ cm^ 

Percent 

error 


Conductivity error of -2.5 percent 


1.1 

26.39 

25.89 

-0.50 

-1.89 

2.1 

29.59 

29.26 

-.33 

-1.11 

3.1 

65.94 

65.30 

-.64 

-.96 

4.1 

170.19 

168.42 

-1.77 

-1.04 

5.1 

276.17 

273.42 

-2.75 

-.99 

6.1 

363.92 

360.58 

-3.34 

-.92 

7.1 

471.32 

467.19 

-4.13 

-.88 

8.1 

595.30 

590.14 

-5.16 

-.87 

9.1 

743.02 

736.54 

-6.48 

-.87 

10.1 

891.44 

883.53 

-7.91 

-.89 

11.1 

1043.89 

1034.39 

-9.50 

-.91 


Conductivity error of -l 

percent 


1.1 

26.39 

25.39 

-1.00 

-3.78 

2.1 

29.59 

28.93 

-.66 

-2.22 

3.1 

65.94 

64.64 

-1.30 

-1.98 

4.1 

170.19 

166.61 

-3.58 

-2.10 

5.1 

276.17 

270.62 

-5.55 

-2.01 

6.1 

363.92 

357,15 

-6.77 

-1.86 

7.1 

471.32 

462.95 

-8.37 

-1.78 

8.1 

595.30 

584.85 

-10.45 

-1.76 

9.1 

743.02 

729.89 

-13.13 

-1.77 

10.1 

891.44 

875.43 

-16.01 

-1.80 

11.1 

1043.89 

1024.71 

-19.18 

-1.84 
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TABLE H.- EFFECT OF ERROR IN CONDUCTIVITY (k) ON 
INVERSE HEAT -FLUX COMPUTATIONS - Concluded 


Time, 

sec 

q 2> 

watts/cm2 

q E’ 

watts/cm^ 

q E " q 2’ 
watts/cm2 

Percent 

error 


Conductivity error of -10 percent 


1.1 

26.39 

24.50 

-1.89 

-7.16 

2.1 

29.59 

29.54 

-.05 

-.17 

3.1 

65.94 

57.50 

-8.44 

-12.80 

4.1 

170.19 

167.24 

-2.94 

-1.73 

5.1 

276.17 

263.98 

-12.19 

-4.41 

6.1 

363.92 

353.82 

-10.10 

- 2.78 

7.1 

471.32 

450.67 

-20.65 

-4.38 

8.1 

595.30 

579.75 

-15.55 

-2.61 

9.1 

743.02 

707.58 

-35.44 

-4.77 

10.1 

891.44 

867.12 

-24.32 

-2.73 

11.1 

1043.89 

998.93 

-44.96 

-4.31 


Conductivity error of +10 percent 


1.1 

26.39 

28.50 

+2.11 

+8.00 

2.1 

29.59 

32.26 

+2.67 

+9.02 

3.1 

65.94 

61.97 

-3.97 

-6.02 

4.1 

170.19 

181.72 

+ 11.53 

+6.77 

5.1 

276.17 

285.72 

+9.55 

+3.46 

6.1 

363.92 

380.74 

+ 16.82 

+4.62 

7.1 

471.32 

483.39 

+12.07 

+2.56 

8.1 

595.30 

621.70 

+26.40 

+4.43 

9.1 

743.02 

759.08 

+16.06 

+2.16 

10.1 

891.44 

930.73 

+39.26 

+4.40 

11.1 

1043.89 

1073.62 

+29.73 

+2.85 



TABLE HI.- EFFECT OF ERROR IN THERMOCOUPLE LOCATION 
ON INVERSE HEAT -FLUX COMPUTATIONS 


Time, 

sec 

<l2’ 

watts/c m2 

^E’ 

watt s/cm 2 

Qe " q 2’ 

watts/cm2 

Percent 

error 

Temperature history at 0.0508 cm assumed to be at 0.0254 cm 

(Error of 0.0254 cm) 

1.1 

26.39 

24.14 

-2.25 

-8.52 

2.1 

29.59 

28.97 

-.62 

-2.07 

3.1 

65.94 

63.26 

-2.68 

-4.06 

4.1 

170.19 

164.46 

-5.73 

-3.37 

5.1 

276.17 

268.73 

-7.44 

-2.69 

6.1 

363.92 

354.78 

-9.14 

-2.51 

7.1 

471.32 

459.70 

-11.62 

-2.47 

8.1 

595.30 

579.29 

-16.01 

-2.69 

9.1 

743.02 

720.41 

-22.61 

-3.04 

10.1 

891.44 

860.61 

-30.83 

-3.46 

11.1 

1043.89 

1001.08 

-42.81 

-4.10 

Temperature history at 0.0762 cm assumed to be at 0.0254 cm 

(Error of 0.0508 cm) 

1.1 

26.39 

21.94 

-4.45 

-16.86 

2.1 

29.59 

28.36 

-1.23 

-4.14 

3.1 

65.94 

60.73 

-5.21 

-7.90 

4.1 

170.19 

158.89 

-11.30 

-6.64 

5.1 

276.17 

261.52 

-14.65 

-5.31 

6.1 

363.92 

346.00 

-17.92 

-4.92 

7.1 

471.32 

448.55 

-22.77 

-4.83 

8.1 

595.30 

564.06 

-31.24 

-5.25 

9.1 

743.02 

699.29 

-43.73 

-5.89 

10.1 

891.44 

832.15 

-59.29 

-6.65 

11.1 

1043.89 

963.44 

-80.45 

-7.71 
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Time, 

sec 


1.1 

2.1 

3.1 

4.1 

5.1 

6.1 

7.1 

8.1 

9.1 

10.1 

11.1 


TABLE IV.- EFFECT OF BIASED TEMPERATURE ON 
INVERSE HEAT -FLUX COMPUTATIONS 


Q2> 

watts/cm2 

<lE’ 

watts/cm^ 

9E ' * 2 ’ 

watts/cm2 

Percent 

error 

Temperature error of +5 percent 


26.39 

27.63 

+1.24 

+4.73 

29.59 

31.41 

+5.34 

+6.14 

65.94 

70.20 

+4.26 

+6.47 

170.19 

180.37 

+10.18 

+5.98 

276.17 

292.29 

+16.12 

+5.84 

363.92 

384.07 

+20.15 

+5.54 

471.32 

495.45 

+24.13 

+5.12 

595.30 

622.55 

+27.25 

+4.58 

743.02 

773.65 

+30.63 

+4.12 

891.44 

923.52 

+32.08 

+3.60 

1043.89 

1078.38 

+34.49 

+3.30 








Temperature, 


>0 


iouu r yuu 


1600 h 800 


> 1400 - £ 700 

•> -p 

H oJ 



cr 

„ 600 


« 500 

d) 


800 L 400 


600 h 300 


Heat flux 
Temperature 


/ 1 / 

/ _ 

ITTffj/ 

• A 

i » / 

t. / cf 

rvi 


/ / 

/ /_ 

/ / / 


/ / / 
/ / / 


“Vi 

v\ H 


'im >k 



i i 

Transition 
100 " begins ^ 


Time, sec 


Figure 4.- Heat-flux history and associated temperature history for different depths in 1.524-cm-thick wall. 
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(a) Input temperature at t = 0.0254 cm. 

Figure 5.- Effect of computing time interval At on computations by the single-thermocouple inverse method, with the exact input for a smoothly 

changing heat flux shown for comparison. = 1.524 cm. 




(b) Input temperature at t = 0.2540 cm. 
Figure 5.- Concluded. 
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(a) At = 0.1 second. 

Figure 6.- Influence of computing interval and temperature-input frequency on oscillations of heat flux computed by the single-thermocouple 
inverse method, with the exact heat-flux input shown for comparison, x = 0.6096 cm: x t = 1.524 cm. 
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(b) At = 0.2 second. 


Figure 6.- Continued. 
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(c) At - 0.4 second. 
Figure 6.- Concluded. 
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Figure 8.- Comparison of single-thermocouple inverse solutions with exact input for an oscillating heat flux, t = 0.0254 cm. 






Heat flux, q, watts/cm' 


Exact neat-flux input 

Single -thermocouple 
inverse solution J 

t = * 025 ^ cm, 

At = .1 sec 


I I 


m K 


i i i 


Transition 

begins 

I i\J 


Jo O O' 



Integral method 

Measurement depths, cm 

T 1 T 2 t 3 Tl i- 

O .025^ .3556 .8636 1.5176 

□ .025^ .5080 1.0160 1.5176 


Time, sec 


Figure 10.- Single-thermocouple inverse solutions and results obtained by the integral method of reference 2 , with the exact heat-flux input 

shown for comparison. = 1.524 cm. 
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Figure 11.- Comparison of exact variation of c^- through wall with integral- method curve fits at 11 sec for q = 1022 watts/cm 2 

and = 1.524 cm. 
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Figure 12.- Solution by method of reference 4 compared with the exact input and with the single-thermocouple inverse solution. 

t = 0.6096 cm ; ij = 1.524 cm. 
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Figure 13.- Single-thermocouple inverse solution obtained by fitting a sixth-order polynomial to scattered temperature data (3a of 25° K) 

compared with exact heat-flux input. 
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Figure 14.- Comparison of heating rates determined from Fire 11 flight temperature measurements during first data period, as computed by single-thermocouple 

inverse method and by integral method in reference 3. 
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(a) = 0.54. 

Figure 16.- Comparison of heating rates determined from Fire II flight temperature measurements during second data period, as computed by 
single-thermocouple inverse method and by integral method in reference 3. = 0.508 cm. 
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Figure 16.- Concluded. 
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